clear all;
clc;

mex cec14_func.cpp;     		% compile the CEC'2014 benchmark set
rng('shuffle');         		% initialization of the random number seed

Dim=50;                	 		% dimension
PS=300;                	 		% population size  
MAX_RUNS=51;            		% epoch
MAX_FES=Dim*10000;      		% maximum fitness evaluations
OUT_TIMES=500;          		% output times of iteration results
u_NS=0.3*PS;            		% upper bound of the neighbor size
c=0.1;                  		% weight value of mean update

function_set=[1:30];    		% index of test functions
function_size = 30;     		% number of test functions
lb=-100;                		% lower bound of test function range
ub=100;                 		% upper bound of test function range

fhd=str2func('cec14_func');  	% load the CEC'2014 benchmark set


for num=1:function_size

    for runs=1:MAX_RUNS

        % record the optimal solution for the first time after the initialization of population and parameters
        pop= lb + (ub-lb)*rand(PS,Dim);             % population initialization
        sigma=0.1;                         			% the fixed variance of random parameters F and CR
        sigma_NS=1;                         		% the fixed variance of random parameters NS
        mu_NS=randi([0.1*PS,0.3*PS]);               % the initial mean value of random parameters NS
        mu_CR=0.5;                          		% the initial mean value of random parameters CR
        mu_F=0.5;                           		% the initial mean value of random parameters F
        A=[];                               		% archieve
        fes=0;                              		% record the number of iterations
        s=1;                                		% Gbest element count
        f= feval(fhd, pop', func_num);      		% individual fitness values
        fes=fes+PS;
        [best_value, index]=min(f);         		% record the optimal solution

        % iteration process
        while fes<MAX_FES
            Sns=[];                         		% save the successful NS values during iteration
            Scr=[];                        		 	% save the successful CR values during iteration
            Sf=[];                          		% save the successful F values during iteration
            NS=Rand_NS(mu_NS,sigma_NS,PS,u_NS);     % random generation of NS based on Cauchy distribution
            F=Rand_F(mu_F,sigma,PS);                % random generation of F based on Cauchy distribution
            CR=Rand_CR(mu_CR,sigma,PS);             % random generation of CR based on Gaussian distribution
            CR=sort(CR);                            % sorting of CR values
            [~,order]=sort(f);                      % sorting of fitness

            for i=1:PS
			
                % match CR and NS for individuals according to the order of individual fitness values
				ns=round(NS(i));
                index_i=find(order==i);
                CR_i=CR(index_i);

                % generate random neighbor region and find the elite individual
                index_neighbor=randperm(PS);     
                optimum_vector=pop(i,:);         
                optimum_index=i;                 
                if ismember(i,index_neighbor(1:ns))
                    index_neighbor(find(index_neighbor(1:ns)==i))=index_neighbor(ns+1);
                end
                [best_value_neighbor,best_index_neighbor]=min(f(index_neighbor(1:ns)));
                if best_value_neighbor<=f(i)
                    optimum_vector=pop(index_neighbor(best_index_neighbor),:);
                    optimum_index=index_neighbor(best_index_neighbor);
                end
                
                % generate the union of archieve and population
                P_A=[pop;A];
                [P_A_PS,~]=size(P_A);

                % select parent individuals for the mutation operation of the target individual
                parent1= unidrnd(PS);
                while(parent1==i||parent1==optimum_index)
                    parent1=unidrnd(PS);
                end
                parent2= unidrnd(P_A_PS);
                while(parent2==i||parent2==parent1||parent2==optimum_index)
                    parent2=unidrnd(P_A_PS);
                end

                % mutation
                v=pop(i,:)+F(i).*(optimum_vector-pop(i,:))+F(i).*(pop(parent1,:)-P_A(parent2,:));

                % cross-border adjustment
                index1 = find(v<lb);
                v(index1)=lb;
                index1 = find(v>ub);
                v(index1)=ub;

                % crossover
                u=crossover_bin(pop(i,:),Dim,CR_i,v);

                % update the fitness of the population
                f1= feval(fhd, u', func_num);
                fes=fes+1;
                if f1<f(i)
                    % keep the archive size at PS
                    [A_line,~]=size(A);
                    if A_line==PS
                        rand_line=unidrnd(PS);                        
                        A(rand_line,:)=pop(i,:);
                    else
                        A=[A;pop(i,:)];
                    end 
                    f(i)=f1;
					% substitution of particles
                    pop(i,:)=u;  
					% save the parameter values of successful individuals
                    Sns=[Sns,NS]; 
                    Scr=[Scr,CR_i];
                    Sf=[Sf,F(i)];
					% save the optimal value
                    if f1<best_value
                        best_value=f1;
                    end
                end
				
                % output optimal solutions
                if fes==s*(MAX_FES/OUT_TIMES)
                  Gbest(s)=best_value;
                  s=s+1;
                end 
            end
			
            % update means
            [~,column]=size(Scr);
            if column~=0
               meanA=sum(Scr)/column;
               meanL=sum(Sf.^2)/sum(Sf);
               meanNS=sum(Sns.^2)/sum(Sns);
            else
               meanA=0;
               meanL=0;
               meanNS=0;
            end
            mu_CR=(1-c)*mu_CR+c*meanA;
            mu_F=(1-c)*mu_F+c*meanL;
            mu_NS=(1-c)*mu_NS+c*meanNS;
        end
    end
end





